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Abstract 

We have calculated evolution of neutron star binaries towards the coalescence 
driven by gravitational radiation. The hydrodynamical effects as well as the general 
relativistic effects are important in the final phase. All corrections up to post 2,5 - 
Newtonian order and the tidal effect are included in the orbital motion. The star is 
approximated by a simple Newtonian stellar model called affine star model. Stellar 
spins and angular momentum are assumed to be aligned. We have showed how the 
internal stellar structure affects the stellar deformation, variations of the spins, and 
the orbital motion of the binary just before the contact. The gravitational wave 
forms from the last a few revolutions significantly depend on the stellar structure. 
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1 Introduction 



Direct detection of gravitational waves will be possible in not so far future by the 
recent constructing laser interferometers such as LIGO [0. One of the promising sources is 
coalescing compact star binary consisting of neutron stars or black holes. The coalescence 
rate of the binary is estimated as ~ 10 yr Mpc -3 [||, |3[]. If we can observe up to a few 
hundred Mpc, we have a few events by a year. The detection of the gravitational wave is 
not only an important probe of general relativity, but also a new tool of measurement of 
Hubble constant and high density nuclear matter 0. ||. 

Theoretically, it is important to understand the final phase of the binary coalescence, 
that is, how the binary evolves from the separation distance of a few hundred km to the 
contact or disruption of the binary. The frequency of the gravitational radiation at this 
phase corresponds to the observational frequency band, about 10 ~ 10 3 Hz. The prepa- 
ration of the gravitational wave form as the template is necessary for the matched filter 
analysis. When the separation L of the binary is much larger than the stellar radius R, 
each star can be regarded as a point mass. The study to get the matched filter is currently 
in progress using the higher order post-Newtonian approximation || [7]]. In particular, the 
test particle limit is extensively studied by the perturbation of the black hole space-time. 
(See e.g., Ref. |S|, |9].) The rotational velocity gradually increases with shrinking of the 
orbital radius, and therefore higher post-Newtonian corrections become crucial. When 
L is a few times R, hydrodynamical effects also become important. The only method 
to examine the final stage may be solving the Einstein equations numerically. Available 
memory and speed of the computers are limited at present, so that the simulations are 
carried out with the Newtonian or post-Newtonian approximation, or preliminary results 
are obtained in fully relativistic case. (See Ref. jlOl, |TT| and therein.) 

One of hydrodynamical effects before the contact of the binary is tidal force. Lai, 
Rasio & Shapiro calculated the dynamical evolution of the binary of simplified stellar 
models under the Newtonian gravity JTJ|, [T3|, |TJj] . They showed that the orbit significantly 
deviates from that of two point masses. The reason is that there is additional contribution 
from the stellar quadrupole moment, which causes marginally stable circular orbit at 
L = (2 ~ 3)R. When the separation becomes below the critical distance, the radial infall 
velocity significantly accelerates in addition to the radiation reaction force. However, 
since the tidal potential depends on L~ 6 , the effect is significant only for close binary, 
i.e., higher post-Newtonian correction. Natural question is whether the contribution is 
still significant, when other post-Newtonian corrections are included. The effect of the 
spinning star, i.e., spin-orbit and spin-spin interaction, may be important since the tidal 
torque changes the spin. They are respectively post L5 -Newtonian and post 2 -Newtonian 
order in the magnitude. The first and second post-Newtonian corrections concerning 
the orbital motion are also important, when the tidal interaction becomes large. In this 
paper, we take all these corrections and the radiation reaction force into consideration, 
and study how the binary evolves at L < 15R. We adopt simplified stellar models 
called affine star models. The model is equivalent to the second order virial equation. 
In this approximation, the fluid displacement is limited to a certain class of the motion, 
i.e., uniform expansion, rotation and quadrupole oscillation. The dynamical degrees of 
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freedom of a star are reduced to finite number, and the equations of motion are not partial 
differential equations, but ordinary differential equations. Using this approximation, we 
can easily simulate the final phase of the binary with various effects concerning spin and 
quadrupole moment. 

In §2, we present our model. We show the collection for the interactions up to the 
post 2 -Newtonian order and the radiation reaction force of post 2,5 -Newtonian order. The 
affine star model is also briefly summarized. We show the deformation of the star, the 
orbital motion, the gravitational wave form and the evolution of the spin near the final 
stage just before the disruption of the binary in §3. In §4, we discuss the implication of 
our results. 



2 Numerical Methods 
2.1 Orbital motion 



In order to solve orbital motion of the binary, we use the Hamiltonian formalism |15 |. 
Regarding the star as a point mass with M a (a = 1, 2), we will examine the relative orbital 
motion. The motion for the center of mass will not be considered here. The equations 
of motion can be described by the total mass Mt and reduced mass /i. We choose the 
orbital plane as the equatorial plane of the polar coordinate (r, 9,4>). The Newtonian, 
the first and second post-Newtonian equations of motion are determined by the following 
Hamiltonians. 

Hp2N = i6cv r + ^J 

GM T J / 2 \ ( 2 P 2 A ,-,,,2,2/ 2 /'••• \ o,,2 i 



3c 4 /z 3 r 



vYr} (3) 



C 2 M 2 \, J 2 pi\ 2 ) G 3 (l + 3v)uM* 

where v = \ijM T . 

Next, we consider the effect of the quadrupole moment of a star. If the star has mass 
quadrupole moment, the Newtonian gravitational potential has an additional contribu- 
tion, ~ 1/r 3 . The relative motion is affected by the monopole-quadrupole interaction 

Ht = "2^ ( 3nV " ^ ( Ma + Ml ih) ^ ' (4) 
where n = r/r and l a is mass quadrupole moment of star a (a = 1, 2). 



3 



The spin-orbit and spin-spin interaction term in appropriate coordinate system can 
be written as fll6|| , 



c 2 r 3 IV 2MJ V 2M : 

^5 = ^3 {3 (Ji • n) (J 2 • n) - (J x • J 2 )} , (6) 

where J a is spin angular momentum of star a and Ljy is orbital angular momentum of 
the Newtonian order, i.e., Ljv = r x p. 

Now, we shall estimate the order of magnitude for these various forces here. Suppose 
that a binary system consists of compact stars with mass M a ~ M and radius R a ~ R 
(a = 1,2). The orbital velocity v at the separation distance L is v ~ [GM / L) 1 / 2 . We 
define the dimensionless relativistic factor, e as e = (GM/Lc 2 ) ~ (v/c) 2 . Since the 
radius of the compact star is of order GM/c 2 , the spin angular momentum J and the 
quadrupole moment I can be written as J = aGM 2 /c and / = (3G 2 M 3 /c 4 , where a and 
/3 are dimensionless numbers. The maximumly rotating Kerr black hole corresponds to 
a = 1. It is likely that a < 1. The value of (3 is small in isolated system, but is induced 
by the companion star as (3 oc (R/L) 3 ~ e 3 . 

We shall write the Newtonian force as /at = GM 2 /L 2 , which can be estimated as 
dHjy/dr ~ Hn/L. The order of the first and second post-Newtonian corrections are ef^ 
and e 2 /jv, respectively. The tidal interaction is the Newtonian force, but the magnitude for 
the compact star is GIM/L 4 ~ /3e 2 fN ~ ^Jn- This is higher post-Newtonian correction. 
The spin-orbit coupling term is formally the post-Newtonian correction, since the term 
contains the factor G/c 2 . The magnitude is, however, estimated as G/c 2 x MLv x J/L 4 ~ 
cte 3 / 2 /jv, that is, the spin-orbit force is the post L5 -Newtonian correction practically. In 
the similar way, the spin-spin interaction is also formally the post-Newtonian correction, 
but practically the post 2 -Newtonian correction as G/c 2 x J 2 /L A ~ a 2 e 2 f^. The equations 
of motion can be determined by the sum of these Hamiltonians. 

2.2 Affine star model 



As seen in the previous section, the orbital motion is affected by the spin angular 
momentum Ji and the mass quadrupole moment i^, as well as the mass M. The spin and 
the quadrupole moment are determined by the stellar structure. We adopted so-called 
the affine star model for the stellar model. In this section, we will summarize the single 
affine star model |T7|, [IS]. The binary consists of two affine star models (a = 1, 2). In this 



section, we omit the suffix a and use the Cartesian coordinate system, whose origin is the 
mass center of a single star. 

The position of a fluid element inside the star, Xi(t) is specified by the Lagrangian 
map from its initial position 

Xi(t) = q ia (t)x a , (7) 
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where q^ a is 3 x 3 matrix and qi a (0) = 5i a . The initial sphere of radius R is transformed 
into an ellipsoid as 

S^XiXj = R\ (8) 

where is the inverse of a symmetric matrix, = qi a qj a - The eigenvalues of 
determine three principal axes of the ellipsoid. 

The dynamics of the fluid motion can be determined in the Hamiltonian formalism 
by the generalized coordinate qi a and its conjugate momentum p ia = I^qia, where the 
inertial moment J is described J = | / x a x a dM and the dot denotes time derivative. 
The initial state is a spherically symmetric star with mass M and radius R, and we 
assume a polytropic equation of state. So that we have J Q = ^K n MR 2 , where K n is a 
constant relating with polytropic index n. The Hamiltonian concerning the internal fluid 
motion is the sum of the kinetic energy T, internal energy U, and self-gravitational energy 
$, 

Hi = T + U + $. (9) 

The kinetic energy is a function of pi a given by 

If 1 

T = 2 J ViVidM = ^fPiaP™- ( 10 ) 

The internal energy is a function of q ia . We use a polytropic equation of state, P = 
Kp 1+( ~ l / n \ then the internal energy of the star is 



U 



f—dM = U \q\- 1/n , (11) 
J o 



nP 
P 

where |q| is determinant of qi a , and Uq is the initial internal energy. The self-gravitational 



energy of the star can be calculated by an elliptical integral |L9 



$ = -$oS tk / '-^du, (12) 

2 Jo \S + u5\* 

where $o is the self-gravitational energy of spherically symmetric star. Since the pressure 
force and self-gravity assumed to balance at the initial state, we have 

$ GM 2 , . 

U = ~ = 77 rp- 13 

3 (5 — n) R 

The quadrupole moment and the spin of the affine star model are given by 

hj = loqiaQja, (14) 

Ji = ^ijkQjaPka- (1^) 
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2.3 Equations of motion and gravitational radiation 



The dynamical degrees of freedom of the system are reduced to finite number in our 
approximation. The variables are fluid ones (qij,Pij) and orbital ones (r, <f>,p r ,P4>)- The 
system is determined by the following total Hamiltonian, 



H = H 



n 



H I2 + H 



N 



H 



PN 



Hp2 



2 N 



H T + H 



so 



H 



ss- 



(16) 



The evolution of 3 a (a = 1, 2) is governed by the familiar formula ||16|| , 



dt 



3GM h 



(n x (I n : n)) 



(17) 



where the subscript b means that 6 = 2 for a = 1, and 6=1 for a = 2. The first term 
of eq.flTT|) represents the precession due to the spin-orbit and spin-spin coupling. This 
term changes the direction, but not the magnitude of the spin. Therefore, when the spins 
and angular momentum are aligned in the same direction, the first term vanishes and the 
direction of the spins is fixed. The second term shows that the magnitude of the spin 
changes through the tidal torque. It is easily checked that the above equation of spin is 
automatically satisfied, if (qij,Pij) follow the equations of motion. Therefore, we do not 
have to calculate the evolution of spins additionally. 

The system (^) is conserved one. We add the effect of the gravitational radiation by 
the lowest quadrupole formula, which is post 2 ' 5 -Newtonian order. The dissipative force 
is introduced by adding the radiation reaction forces (f r , f^, F r , F^) in the Hamiltonian 
formalism (e.g., Ref.[^, |21|). The explicit radiation reaction forms are given by 



fr 



Fr 

Fa, 



8G 2 



15c 5 ur 2 

8G 2 PrPj, 
3c 5 z/r 4 
8G 2 p r f 

' 3c 5 r A 



2p r 2 + 



GuM 3 " 



vr 



SG 2 P4> f 2Gv 2 M 3 | 2p\ 2 N 
5c 5 vr 3 \ r r 2 ^ r , 



(18) 

(19) 
(20) 

(21) 



The evolution of the variables is determined by the following set of equations, 

dH 



dr 
~dt 
d<\> 
~dt 

dp r 

~dt 



dp r 

dH 

dp^ 
dH 
dr 



+ fr 



+ 



+ F r 



(22) 
(23) 
(24) 
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dt 
dqij 

dt 
dpij 

dt 



dH 
~ 30 
dH 
dp^ ' 
dH 
dqij ' 



+ F 6 



(25) 
(26) 
(27) 



We describe two polarization modes of gravitational waves by h + and h x . They are 
written by the distance to the source D and the directional angle 0. Including the 
quadrupole moment terms of each ellipsoids, the polarization modes are given by, 



G 



2c 4 D 

dt 2 
G 



(i + cos 2 e) 



2c 4 D 
cP_ 
dt 2 



fir 2 cos 20+ £ {(Ia)i-(Ia) 2 }cos2(0 + ^ 

a=l,2 

(2 cos 6) 

/ir 2 sin20 + {(Ia)i- (Ia)2}sin2(0 + a 

a=l,2 



(2f 



(29) 



where (I„ 



1,2) is the inertial moment of the ellipsoid along the principal axis ctj 



Both ai and a,2 axes are chosen on the orbital plane of the binary. The a± axis is initially 
chosen as the separation direction n of the binary and the a2 axis is orthogonal to it. 
Three eigen- values of the matrix (I a )ij give (I a )i {i = 1,2,3). In eqs.(|28|), (p9"D, 9 a is the 
lag angle, which is defined by the angle between the major axis of the ellipsoid a and the 
direction to the companion star. Tidally locked state corresponds to 6 a = 0. 



3 Results 



We have calculated the evolution of close binary from a circular orbit of the separa- 
tion r = 15R to the contact, using two identical affine stars. The most natural orbit 
near r = 15R is circular one, because the eccentricity is diminished by the gravitational 
radiation as the binary orbit shrinks. As for the internal structure at r = 15R, we 
adopt the equilibrium solution of the ellipsoid, i.e., irrotational Roche- Riemann ellipsoid 
or Maclaurin- Jacobi- Roche ellipsoid [^2| . The former solution basically represents tidally 
locked, static figure in the inertial frame. We adopt Maclaurin- Jacobi-Roche solution as 
the spinning stellar model. Both stellar spins are assumed to be parallel to the orbital 
angular momentum, but the initial magnitude of spin is variable. The magnitude is pa- 
rameterized by x — |Jo| /(I<j)3^o at the initial separation, where fl Q is the orbital angular 
frequency of the star. We selected the parameter as x — —2,-1,-0.5,0,0.5,1,2. The 
maximum rotation \x\ = 2 in our model corresponds to the period of ~ AOnto, where to 
is the dynamical time at an individual star surface i.e., tr, = (i? 3 /GM) 1//2 . Typically, it 
is about ten milliseconds for the actual neutron star. The results given below are almost 
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for x — except in §3.3 and §3.4, because the results for \ 7^ models are almost the 
same as those for x — 0- The stellar model is set up as (GM/Rc 2 ) ~ 0.2 with polytropic 
index n = 0.5 or n — 1, corresponding to the typical neutron star of M ~ 1.4M and 
i? ~ 10 km. We have numerically calculated the subsequent evolution from these initial 
conditions. We ignored the dissipative viscosity terms for internal motions during the 
evolution of the orbit, so that the stellar circulation is conserved before tidal torque is 
important. We will give our results of numerical simulations below. 

3.1 Tidal distortion 

We show the typical examples of the quadrupole deformation of stars with n = 0.5 
and n — 1 in Fig.l. The initial state corresponds to x = 0, but the results are almost the 
same for x 7^ 0. Three orthogonal principal axes, a±, a<i and 03, of the ellipsoid are shown 
as a function of the orbital separation. The minimum separation distance corresponds to 
the contact of the binary. That is, r ~ 2AR for n = 0.5 and r ~ 2.3R for n — 1. The 
axes ai and a 2 are located on the orbital plane, while a 3 is the direction of the orbital 
angular momentum. The ai-axis is initially chosen as the separation direction n of the 
binary, but the direction does not exactly follow the orbital motion during the dynamical 
evolution. The orbital revolution at small r becomes so fast that the principal axis slightly 
deviates from the direction n, i.e., there is a time lag. We will discuss the lag angle in 
next subsection. However, the lag angle is small, so that the star is elongated almost 
in the separation direction n and compressed in other two directions. The quadrupole 
deformation becomes significant for r < 5R and grows as \ai/R — 1| oc r~ 3 . The internal 
structure also affects the shape of the ellipsoids. The model with n — 1 is more centrally 
condensed than that with n = 0.5, so that the shape is less deformed. The difference is, 
however, a few percent in the magnitude. We will see how this difference appears in the 
orbital motion and gravitational radiation below. 

3.2 Orbital motions 

The binary evolution driven by the gravitational radiation is calculated in terms of four 
approximations to estimate various forces. (1) All post 2 -Newtonian, spin and tidal effects 
are included. (2) The post 2 -Newtonian effects are included, but the stars are regarded as 
point masses. (3) Purely Newtonian gravity is used, but the tidal effects are included. 
(4) The Newtonian gravity is used, and the stars are regarded as point masses. We show 
the angular and radial velocities in Fig. 2 and Fig.3. We give the comparison between 
the approximations (1) and (2) in Fig. 2a and Fig.3a, and the comparison between the 
approximations (3) and (4) in Fig.2b and Fig. 3b. The difference due to the approximations 
becomes clear below r ~ 3R, where the star is much deformed. The model with n = 0.5 is 
more deformed and therefore the deviation from the point mass approximation becomes 
more significant. Although the spin-orbit force with positive x are repulsive, the post- 
Newtonian corrections, spin-spin and tidal force are dominant attractive forces. Combined 
with all these effects, the final approaching velocity significantly accelerates. 
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The absolute magnitude of the velocity shifts, when the post-Newtonian corrections 
are included. The ratio, i.e., the direction of colliding stars, may be more useful than 
the absolute values as the initial conditions of the hydrodynamical simulation after the 
coalescence. We show the velocity ratio —V r /V,p in Fig. 4. The approaching velocity 
increases up to ten percent of the orbital angular velocity, i.e., the velocity ratio — V r /V^ ~ 
0.1 at the contact of the binary. This ratio further increases with decrease of the polytropic 
index n, because the tidal force is larger. 

When the binary separation is large enough, the binary evolves at much longer time- 
scale compared with that of the internal motions. The stellar figure is tidally locked 
at this stage. As the time-scale of the orbital motion decreases, the tidal lock becomes 
loose and the lag angle becomes non-zero. We show the lag angles for two stellar models 
(n = 1 and n = 0.5) in Fig.5. The negative lag angles means that the elliptical figure 
does not follow the orbital revolution, and that there is a time lag to settle down. The 
lag angles amount to 0.1 radian just before the contact. This suggests that both stars at 
the coalescence touch each other not in the longest axis of the ellipsoid, but the direction 
slightly deviates. 

3.3 Gravitational wave forms 

The tidal force affects the orbital motion just before the contact. We have synchronized 
the phase of the gravitational wave at the initial separation r = 15-R. In Fig. 6 and Fig. 7, 
we only show the evolution of the gravitational wave forms radiated at the last several 
revolutions. The correction of tidal force in the amplitude is of order I / (/xr 2 ) ~ f3e 2 . 
Therefore, tidal effect is a small correction in the amplitude, but it accumulates in the 
phase of the wave. As the time proceeds, the wave form deviates from that of point 
masses, which is represented by the dotted line in Fig. 6. The deviation of the n = 0.5 
model is larger than that of the n = 1 model, because of larger tidal deformation. The 
phase shift becomes significant only for the last a few revolutions, and the amount of the 
phase shift is less than 7r, if there is no initial spin. In Fig. 7, we show the dependence of 
the initial spin. The stellar spin gives large difference in the phase shift. For the larger 
initial spin, the repulsive spin-orbit force becomes stronger, so that it takes much longer 
time to contact each other. The amount of the phase lag becomes more than 7r for the 
fast spinning case. 

3.4 Evolution of spin 

The stellar spin is modified through the tidal torque following eq.([T7|). We have 
calculated the evolution of the spin from various initial states. We show the typical 
evolution of the spin in Fig.8. The variations of spins are quite small at large separation, 
but monotonically grow only in the final stage, r < 4R. The stars significantly deviate 
from the spherical state as mentioned in §3.1, and the tidal torque becomes quite large 
there. Therefore, variation of the spin is quite a steep function of r, since the tidal torque 
is dJ/dt oc Ir~ 3 oc r~ 6 . We also show the total amount of the change of the spin for the 
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various initial spin angular momentum in Fig.9. The transferee! momentum through the 
tidal torque decreases with the increase of the parameter \. The angular momentum is 
much more transfered in the counter-rotating star with the orbital angular momentum, 
i.e, star with negative x- 



4 Discussion 



The simple estimate of the tidal effect is of order ~ P^Jn ~ £ 5 I'n, as shown in 
§2.1. The final state of the binary corresponds to e ~ 0.3, so that the correction seems 
to be (0.3) 5 ~ 0.002, but is much large in the actual problem. The discrepancy comes 
from the existence of the marginal circular radius, where the hydrodynamical instability 
occurs [|T2 . If the radiation reaction force is not included, there are circular orbits for any 



separations, but they are not always stable |23]. We expect that there exists innermost 
stable circular orbit, which the post-Newtonian gravity and the tidal force cause at the 
critical separation. Two stars significantly accelerate in the radial direction below this 
critical radius. In a future paper, we will show that the location of the innermost stable 
circular orbit is affected by polytropic index n and initial spin parameter x- The tidal 
effects at the last a few cycles cause large radial infall velocity, lag angle and spin-up. 
The higher post-Newtonian corrections or fully relativistic hydrodynamical calculation 
will be required to study much more correct final phase. The topic treated in this paper 
is preliminary, but will be useful as the guide for the future theoretical work and the 
observations by the advanced type of the laser interferometer. 
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Figure Captions 



Fig. la: Variations of the tri-axes of ellipsoid with n = 0.5 and x — as a function of 
the binary separation. The solid, dashed, and dotted lines indicate the axis a±, a 2 , and a 3 . 

Fig. lb: The same as Fig. la, but with n — 1. 

Fig. 2a: The angular velocity at the binary separation r. The solid line denotes the 
result included post 2 -Newtonian, spin-orbit, spin-spin and tidal forces with n = 0.5 stel- 
lar model without initial spin. The dashed line is the same but with n = 1 model. For 
comparison, the result of point mass, neglecting the tidal and spin effects, is described by 
the dotted line. 

Fig. 2b: The same as Fig. 2a, but the results are based on the purely Newtonian cal- 
culation. The solid and dashed line denote the results included tidal forces with n = 0.5 
and n = 1 models, respectively. The dotted line is the result of point mass. 

Fig. 3a: The radial velocity of the binary with relativistic corrections. The meaning 
of each line is the same as in Fig. 2a. 

Fig. 3b: The radial velocity of the binary with Newtonian gravity. The meaning of 
each line is the same as in Fig. 2b. 

Fig. 4a: The velocity ratio — V r jV^ in the coalescing binary with relativistic corrections. 
The meaning of each line is the same as in Fig. 2a. 

Fig. 4b: The velocity ratio — V r /V^, in the coalescing binary with Newtonian gravity. 
The meaning of each line is the same as in Fig. 2b. 

Fig. 5: The lag angle for the binary without the initial spin. The solid and dashed 
lines correspond to n — 0.5 and n — 1, respectively. 

Fig. 6a: The gravitational wave forms h + from the binary with polytropic index n. The 
solid and dashed lines correspond to the result of n — 0.5 and that of n = 1, respectively. 
Result of point masses is also given by the dotted line for comparison. 

Fig. 6b: The same as Fig. 6a, but for the gravitational wave forms h x . 

Fig. 7a: The gravitational wave forms h + of the binary with n = 0.5 from various initial 
spin states. The solid, dotted and dashed lines correspond to x = 0, 1, 2. 

Fig. 7b: The same as Fig.7a, but the initial spin is counter-rotating. The solid, dot- 
ted and dashed lines correspond to x = 0, — 1, —2. 

Fig. 8: Variation of stellar spin with all corrections. The solid line represents the evolu- 
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tion for n = 0.5 model, and the dashed line for n = 1. The initial spin is zero. 

Fig. 9: Total change of the spin angular momentum for various initial x values. The 
filled circles and triangles correspond to n — 0.5 and n — 1 model. 
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